Wprowadzenie
Na przestrzeni ostatnich lat zauważono stopniowy spadek rozmiaru śledzia oceanicznego wyławianego w Europie. Do analizy zebrano pomiary śledzi i warunków w jakich żyją z ostatnich 60 lat. Dane były pobierane z połowów komercyjnych jednostek. W ramach połowu jednej jednostki losowo wybierano od 50 do 100 sztuk trzyletnich śledzi.
Inicjalizacja środowiska
library(ggplot2)
library(dplyr)
library(DT)
library(GGally)
library(magick)
library(cowplot)
# ustawmy seed'a dla randoma w celu powtarzalności wyników
set.seed(1234)
Zainicjujmy również funkcje pomocnicze do prezentacji danych
prettyTable <- function(table_df, round_columns=numeric(), round_digits=2) {
DT::datatable(table_df, style="bootstrap", filter = "top", rownames = FALSE, extensions = "Buttons",
options = list(dom = 'Bfrtip', buttons = c('copy', 'csv', 'excel', 'pdf', 'print'))) %>%
formatRound(round_columns, round_digits)
}
Import danych
Dane znajdują się w pliku sledzie.csv załączonym w repozytorium. Opisy ich atrybutów można znaleźć w treści zadania; dla wygody jednak zamieszczono go również poniżej:
Kolejne kolumny w zbiorze danych to:
- length: długość złowionego śledzia [cm];
- cfin1: dostępność planktonu [zagęszczenie Calanus finmarchicus gat. 1];
- cfin2: dostępność planktonu [zagęszczenie Calanus finmarchicus gat. 2];
- chel1: dostępność planktonu [zagęszczenie Calanus helgolandicus gat. 1];
- chel2: dostępność planktonu [zagęszczenie Calanus helgolandicus gat. 2];
- lcop1: dostępność planktonu [zagęszczenie widłonogów gat. 1];
- lcop2: dostępność planktonu [zagęszczenie widłonogów gat. 2];
- fbar: natężenie połowów w regionie [ułamek pozostawionego narybku];
- recr: roczny narybek [liczba śledzi];
- cumf: łączne roczne natężenie połowów w regionie [ułamek pozostawionego narybku];
- totaln: łączna liczba ryb złowionych w ramach połowu [liczba śledzi];
- sst: temperatura przy powierzchni wody [°C];
- sal: poziom zasolenia wody [Knudsen ppt];
- xmonth: miesiąc połowu [numer miesiąca];
- nao: oscylacja północnoatlantycka [mb].
Wiersze w zbiorze są uporządkowane chronologicznie.
sourceDataFileName <- "sledzie.csv"
Dokonajmy wstępnego wczytania danych - pierwszych 100 wierszy w celu automatycznego określenia klasy oraz podglądu danych
initialHerringData <- read.csv(sourceDataFileName, nrows = 100)
classes <- sapply(initialHerringData, class)
classes
X length cfin1 cfin2 chel1 chel2 lcop1 lcop2 fbar recr cumf totaln sst sal
"integer" "numeric" "factor" "factor" "factor" "factor" "factor" "factor" "numeric" "integer" "numeric" "numeric" "factor" "numeric"
xmonth nao
"integer" "numeric"
Jak widać, większość danych została zakwalifikowana jako factor mimo, że z opisu wynika, że są to liczby. Wczytajmy zatem cały plik w odpowiednich klasach - wartości brakujące reprezentuje ?, natomiast jako separator dziesiętny przyjmujemy kropkę .. Plik nie zawiera również komentarzy. Szybka analiza pliku pozwoliła również na identyfikację, że zbiór składa się z około 52.6 tys wierszy - informację tą wykorzystamy w celu przyśpieszenia importu.
classes <- c("NULL", rep(c("numeric"), 13), "factor", "numeric")
herringData <- read.csv(sourceDataFileName, comment.char = "", header=TRUE, colClasses = classes, nrows=52600, dec = ".", na.strings = c("?"))
Omówienie danych
Sprawdzmy czy rozmiary zbioru pokrywają się z oczekiwaniami
dim(herringData)
[1] 52582 15
Wygląda na to, że dane zostały poprawnie wczytane. Zobaczmy je więc
prettyTable(head(herringData, n = 100))
Następnie przeprowadzmy krótkie podsumowanie
summary(herringData)
length cfin1 cfin2 chel1 chel2 lcop1 lcop2 fbar
Min. :19.0 Min. : 0.0000 Min. : 0.0000 Min. : 0.000 Min. : 5.238 Min. : 0.3074 Min. : 7.849 Min. :0.0680
1st Qu.:24.0 1st Qu.: 0.0000 1st Qu.: 0.2778 1st Qu.: 2.469 1st Qu.:13.427 1st Qu.: 2.5479 1st Qu.:17.808 1st Qu.:0.2270
Median :25.5 Median : 0.1111 Median : 0.7012 Median : 5.750 Median :21.673 Median : 7.0000 Median :24.859 Median :0.3320
Mean :25.3 Mean : 0.4458 Mean : 2.0248 Mean :10.006 Mean :21.221 Mean : 12.8108 Mean :28.419 Mean :0.3304
3rd Qu.:26.5 3rd Qu.: 0.3333 3rd Qu.: 1.7936 3rd Qu.:11.500 3rd Qu.:27.193 3rd Qu.: 21.2315 3rd Qu.:37.232 3rd Qu.:0.4560
Max. :32.5 Max. :37.6667 Max. :19.3958 Max. :75.000 Max. :57.706 Max. :115.5833 Max. :68.736 Max. :0.8490
NA's :1581 NA's :1536 NA's :1555 NA's :1556 NA's :1653 NA's :1591
recr cumf totaln sst sal xmonth nao
Min. : 140515 Min. :0.06833 Min. : 144137 Min. :12.77 Min. :35.40 8 : 9920 Min. :-4.89000
1st Qu.: 360061 1st Qu.:0.14809 1st Qu.: 306068 1st Qu.:13.60 1st Qu.:35.51 10 : 7972 1st Qu.:-1.89000
Median : 421391 Median :0.23191 Median : 539558 Median :13.86 Median :35.51 7 : 6922 Median : 0.20000
Mean : 520367 Mean :0.22981 Mean : 514973 Mean :13.87 Mean :35.51 9 : 5714 Mean :-0.09236
3rd Qu.: 724151 3rd Qu.:0.29803 3rd Qu.: 730351 3rd Qu.:14.16 3rd Qu.:35.52 6 : 4218 3rd Qu.: 1.63000
Max. :1565890 Max. :0.39801 Max. :1015595 Max. :14.73 Max. :35.61 5 : 3736 Max. : 5.08000
NA's :1584 (Other):14100
Jak widać, śledzie mają długość w przedziale 19-32.5 cm, z medianą ustaloną dla wartości 25.5 cm, najwięcej śledzi wyławia się w sierpniu, październiku oraz lipcu. Można również stwierdzić, że liczba brakujących danych dla pól cfin, chel1, lcop oraz sst jest porównywalna; dla innych column wszystkie dane są określone. Możemy również zwrócić uwagę na wielkość połowów - ich liczba mieści się w przedziale 14.4 - 101.6 tys z medianą na poziomie 54 tys oraz średnią 51.5 tys. Rozkład ten jest więc prawdopodobnie zbliżony do normalnego, co z resztą możemy sprawdzić.
ggplot(data=herringData, aes(x=totaln)) +
geom_density() +
labs(x = "Wielkość połowu", y = "Gęstość") +
theme_bw()

Niekoniecznie - nie przypomina on żadnego z popularnych rozkładów.
Mamy za zadanie zweryfikować, czy przez ostatnie kilka lat rozmiar śledzia spadł. Przydatny może okazać sie więc wykres prezentujący średni rozmiar na rok. Jest jednak jeden problem - w danych mamy jedynie miesiąc połowu, informację o tym, że są to dane chronologicznie posortowane oraz pochodzą z ostatnich 60 lat. Należy więc w miarę możliwości dodać rok do danych bazując na tych informacjach. Sprawdzmy jednak czy nie ma wśród tych danych brakujących wartości
sum(is.na(herringData$xmonth))
[1] 0
brak - możemy więc przejść dalej. Dane są chronologiczne, także powinniśmy być w stanie wywnioskować rok na podstawie zmiany miesiąca na mniejszy
evaluateYearByMonthChange <- function(fishing) {
curYear <- 0
recentMonth <- 0
l <- lapply(fishing, function(v) {
v<- as.integer(v)
if (v < recentMonth)
curYear <<- curYear + 1
recentMonth <<- v
curYear
})
unlist(l)
}
herringDataWithYearByMonthChange <- mutate(herringData, year=evaluateYearByMonthChange(xmonth))
herringDataWithYearByMonthChange
Pogrupujmy je po roku dzięki funkcji
statsInYear <- function(dataWithYear) {
dataWithYear %>%
group_by(year) %>%
summarise_at(vars(length), funs(mean(., na.rm=TRUE), min = mean - sd(., na.rm= TRUE), max = mean + sd(., na.rm= TRUE)))
}
meanLengthInYearChart <- function(dataWithYear) {
meanPerYear <- statsInYear(dataWithYear)
ggplot(data=meanPerYear, aes(x = year, y = mean)) +
geom_ribbon(aes(ymin = min, ymax = max), alpha = 0.5, fill = "darkseagreen3", color = "transparent") +
geom_line(color = "aquamarine4", lwd = 0.7) +
labs(x = "Rok", y = "Długość śledzia") +
theme_bw()
}
meanLengthInYearChart(herringDataWithYearByMonthChange)

Założenie o chronologii jak widać okazało się błędne, przynajmniej w zakresie kolejności miesięcy. Jest również pole recr mówiące o połowie w roku. Widzimy, ze występuje pewna powtarzalność, więc jeśli dane byłyby chronologiczne, a zarazem ich rozpiętość jest całkiem spora, powinno to pozwolić na zgrupowanie.
evaluateYear <- function(fishing) {
curYear <- 0
recentYearFishing <- 0
l <- lapply(fishing, function(v) {
if (v != recentYearFishing)
curYear <<- curYear + 1
recentYearFishing <<- v
curYear
})
unlist(l)
}
herringDataWithYearByRecr <- mutate(herringData, year=evaluateYear(recr))
herringDataWithYearByRecr
meanLengthInYearChart(herringDataWithYearByRecr)

Niestety, ponownie porażka. W takim razie może po prostu sprawdzmy czy występuje trend w czystych danych:
ggplot(herringData, aes(x=as.numeric(row.names(herringData)), y=herringData$length)) +
geom_line(color = "aquamarine4") +
labs(x = "Rok", y = "Długość śledzia") +
theme_bw()

Dobrze, jakąś zależność widać - może aby być jej pewnym uśrednijmy dane stosując prawo wielkich liczb oraz wiedzę o tym, że mamy do czynienia z danymi z okresu 60 lat
expectedYears <- 60
expectedRowsPerYears <- nrow(herringData) %/% 60
herringDataWithYearByPeriod <- herringData %>%
mutate(year = 1:n() %/% expectedRowsPerYears)
herringDataWithYearByPeriod
meanLengthInYearChart(herringDataWithYearByPeriod)
Error in mean.default(length, na.rm = TRUE, function (x, na.rm = FALSE, :
'trim' must be numeric of length one
Widzimy więc, że w rzeczywistości istnieje istotny spadek średniej długości śledzia na przełomie lat. Możemy to zwizualizować.
img <- image_read("herring.svg")
img <- image_trim(img)
heightRatio <- 110/475
maxWidth = 542
maxHeight <- heightRatio * maxWidth
herringsImgs <- statsInYear(herringDataWithYearByPeriod) %>%
apply(1, function(r) {
newWidth <- round((r[[2]] - 10) * 30)
widthDif <- round((maxWidth - newWidth)/2)
heightDif <- round(widthDif * heightRatio)
image_scale(img, newWidth) %>%
image_border("white", paste(c(widthDif, "x", heightDif), collapse = "")) %>%
image_annotate(paste(c("Rok: ", r[[1]]), collapse = ""), color = "black", size = 30)
})
frames <- image_morph(image_join(herringsImgs), frames = 1)
anim <- image_animate(frames, fps = 10)
image_write(anim, "herring.gif")

Zidentyfikujmy możliwe powody takiego zachowania wyznaczając macierz korelacji.
Wcześniej jednak należy zmienić typ kolumny xmonth z powrotem na liczbę i przygotować dane.
herringDataWithNumericMonth = mutate(herringData, xmonth = as.numeric(xmonth))
Zobaczmy wynik w bardziej przyjaznej formie
ggcorr(
herringDataWithNumericMonth,
nbreaks = 9,
label = TRUE,
label_size = 3,
color = "grey50",
layout.exp = 1,
hjust = 0.75,
)

Widać dużą korelację pomiędzy parami chel1 - lcop1, chel2 - lcop2, fbar - cumf, cfin2 - lcop2 oraz cumf - totaln. Na samą długość w największym stopniu wpływa sst, nao oraz fbar. Widzimy że lcop1, lcop2 oraz cumf mogą zostać pominięte w procesie dalszej analizy jako redundantne. Teoretycznie moglibyśmy również usunąć miesiąc ze względu na powtarzalność, ale może mieć on wpływ na cykl rozwojowy, więc pozostawmy go póki co.
filteredHerringDataWithNumericMonth <- select(herringDataWithNumericMonth, -c(lcop1, lcop2, cumf))
Zobaczmy jeszcze jak przebiega zależność pomiędzy najbardziej obiecującymi zmiennymi na długość; kolorami wyróżniono również miesiące
numericChart <- function(colName) {
ggplot(herringData, aes(x = length, y = herringData[[colName]], color=xmonth)) +
geom_jitter(size = 1.5, stat = "identity") +
labs(x = "Długość", y = colName)
}
plot_grid(
numericChart("chel1"),
numericChart("fbar"),
numericChart("sst"),
numericChart("nao")
)
Powiązanie rosnącego sst oraz nao jasno wpływa negatywnie na długość, przy fbar również można dostrzec pozywytną korelację.
Sprawdzmy również zależność pomiędzy miesiącem połowu a długością; w macierzy korelacji nie mogliśmy tego do konca zweryfikować, jednak ze względu na cykliczność było to również utrudnione (miesiąc w tym wypadku jest sztucznie liczbą).
numericChart("xmonth")

Wydaje się, że nie ma żadnej znaczącej zależności pomiędzy miesiącem połowu a długością śledzia
LS0tDQp0aXRsZTogIlJhcG9ydCB6IGFuYWxpenkgem1pYW55IHcgcm96bWlhcnplIHd5xYJhd2lhbmVnbyDFmmxlZHppYSBPY2Vhbmljem5lZ28iDQphdXRob3I6ICJXaXRvbGQgS3VwxZsiDQpvdXRwdXQ6IA0KICBodG1sX25vdGVib29rOg0KICAgIHRvYzogdHJ1ZQ0KICAgIHRvY19mbG9hdDogdHJ1ZQ0KICBnaXRodWJfZG9jdW1lbnQ6IGRlZmF1bHQNCmRhdGU6ICJgciBmb3JtYXQoU3lzLnRpbWUoKSwgJyVkICVCLCAlWScpYCINCi0tLQ0KIyBBYnN0cmFrdA0KDQojIFdwcm93YWR6ZW5pZQ0KDQpOYSBwcnplc3RyemVuaSBvc3RhdG5pY2ggbGF0IHphdXdhxbxvbm8gc3RvcG5pb3d5IHNwYWRlayByb3ptaWFydSDFm2xlZHppYSBvY2Vhbmljem5lZ28gd3nFgmF3aWFuZWdvIHcgRXVyb3BpZS4gRG8gYW5hbGl6eSB6ZWJyYW5vIHBvbWlhcnkgxZtsZWR6aSBpIHdhcnVua8OzdyB3IGpha2ljaCDFvHlqxIUgeiBvc3RhdG5pY2ggNjAgbGF0LiBEYW5lIGJ5xYJ5IHBvYmllcmFuZSB6IHBvxYJvd8OzdyBrb21lcmN5am55Y2ggamVkbm9zdGVrLiBXIHJhbWFjaCBwb8WCb3d1IGplZG5laiBqZWRub3N0a2kgbG9zb3dvIHd5YmllcmFubyBvZCA1MCBkbyAxMDAgc3p0dWsgdHJ6eWxldG5pY2ggxZtsZWR6aS4NCg0KIyMjIEluaWNqYWxpemFjamEgxZtyb2Rvd2lza2ENCmBgYHtyLCBlY2hvID0gVFJVRSwgcmVzdWx0cz0naGlkZSd9DQpsaWJyYXJ5KGdncGxvdDIpDQpsaWJyYXJ5KGRwbHlyKQ0KbGlicmFyeShEVCkNCmxpYnJhcnkoR0dhbGx5KQ0KbGlicmFyeShtYWdpY2spDQpsaWJyYXJ5KGNvd3Bsb3QpDQojIHVzdGF3bXkgc2VlZCdhIGRsYSByYW5kb21hIHcgY2VsdSBwb3d0YXJ6YWxub8WbY2kgd3luaWvDs3cNCnNldC5zZWVkKDEyMzQpDQpgYGANCg0KWmFpbmljanVqbXkgcsOzd25pZcW8IGZ1bmtjamUgcG9tb2NuaWN6ZSBkbyBwcmV6ZW50YWNqaSBkYW55Y2gNCmBgYHtyfQ0KcHJldHR5VGFibGUgPC0gZnVuY3Rpb24odGFibGVfZGYsIHJvdW5kX2NvbHVtbnM9bnVtZXJpYygpLCByb3VuZF9kaWdpdHM9Mikgew0KICAgIERUOjpkYXRhdGFibGUodGFibGVfZGYsIHN0eWxlPSJib290c3RyYXAiLCBmaWx0ZXIgPSAidG9wIiwgcm93bmFtZXMgPSBGQUxTRSwgZXh0ZW5zaW9ucyA9ICJCdXR0b25zIiwNCiAgICAgICAgICAgICAgICAgIG9wdGlvbnMgPSBsaXN0KGRvbSA9ICdCZnJ0aXAnLCBidXR0b25zID0gYygnY29weScsICdjc3YnLCAnZXhjZWwnLCAncGRmJywgJ3ByaW50JykpKSAlPiUNCiAgICBmb3JtYXRSb3VuZChyb3VuZF9jb2x1bW5zLCByb3VuZF9kaWdpdHMpDQp9DQpgYGANCg0KIyMjIEltcG9ydCBkYW55Y2gNCg0KRGFuZSB6bmFqZHVqxIUgc2nEmSB3IHBsaWt1IGBzbGVkemllLmNzdmAgemHFgsSFY3pvbnltIHcgcmVwb3p5dG9yaXVtLiBPcGlzeSBpY2ggYXRyeWJ1dMOzdyBtb8W8bmEgem5hbGXFusSHIHcgW3RyZcWbY2kgemFkYW5pYV0oaHR0cDovL3d3dy5jcy5wdXQucG96bmFuLnBsL2FsYWJpamFrL2VtZC9wcm9qZWt0L3Byb2pla3RfYW5hbGl6YS5odG1sKTsgZGxhIHd5Z29keSBqZWRuYWsgemFtaWVzemN6b25vIGdvIHLDs3duaWXFvCBwb25pxbxlajoNCg0KS29sZWpuZSBrb2x1bW55IHcgemJpb3J6ZSBkYW55Y2ggdG86DQoNCi0gbGVuZ3RoOiBkxYJ1Z2/Fm8SHIHrFgm93aW9uZWdvIMWbbGVkemlhIFtjbV07DQotIGNmaW4xOiBkb3N0xJlwbm/Fm8SHIHBsYW5rdG9udSBbemFnxJlzemN6ZW5pZSBDYWxhbnVzIGZpbm1hcmNoaWN1cyBnYXQuIDFdOw0KLSBjZmluMjogZG9zdMSZcG5vxZvEhyBwbGFua3RvbnUgW3phZ8SZc3pjemVuaWUgQ2FsYW51cyBmaW5tYXJjaGljdXMgZ2F0LiAyXTsNCi0gY2hlbDE6IGRvc3TEmXBub8WbxIcgcGxhbmt0b251IFt6YWfEmXN6Y3plbmllIENhbGFudXMgaGVsZ29sYW5kaWN1cyBnYXQuIDFdOw0KLSBjaGVsMjogZG9zdMSZcG5vxZvEhyBwbGFua3RvbnUgW3phZ8SZc3pjemVuaWUgQ2FsYW51cyBoZWxnb2xhbmRpY3VzIGdhdC4gMl07DQotIGxjb3AxOiBkb3N0xJlwbm/Fm8SHIHBsYW5rdG9udSBbemFnxJlzemN6ZW5pZSB3aWTFgm9ub2fDs3cgZ2F0LiAxXTsNCi0gbGNvcDI6IGRvc3TEmXBub8WbxIcgcGxhbmt0b251IFt6YWfEmXN6Y3plbmllIHdpZMWCb25vZ8OzdyBnYXQuIDJdOw0KLSBmYmFyOiBuYXTEmcW8ZW5pZSBwb8WCb3fDs3cgdyByZWdpb25pZSBbdcWCYW1layBwb3pvc3Rhd2lvbmVnbyBuYXJ5Ymt1XTsNCi0gcmVjcjogcm9jem55IG5hcnliZWsgW2xpY3piYSDFm2xlZHppXTsNCi0gY3VtZjogxYLEhWN6bmUgcm9jem5lIG5hdMSZxbxlbmllIHBvxYJvd8OzdyB3IHJlZ2lvbmllIFt1xYJhbWVrIHBvem9zdGF3aW9uZWdvIG5hcnlia3VdOw0KLSB0b3RhbG46IMWCxIVjem5hIGxpY3piYSByeWIgesWCb3dpb255Y2ggdyByYW1hY2ggcG/Fgm93dSBbbGljemJhIMWbbGVkemldOw0KLSBzc3Q6IHRlbXBlcmF0dXJhIHByenkgcG93aWVyemNobmkgd29keSBbwrBDXTsNCi0gc2FsOiBwb3ppb20gemFzb2xlbmlhIHdvZHkgW0tudWRzZW4gcHB0XTsNCi0geG1vbnRoOiBtaWVzacSFYyBwb8WCb3d1IFtudW1lciBtaWVzacSFY2FdOw0KLSBuYW86IG9zY3lsYWNqYSBww7PFgm5vY25vYXRsYW50eWNrYSBbbWJdLg0KDQpXaWVyc3plIHcgemJpb3J6ZSBzxIUgdXBvcnrEhWRrb3dhbmUgY2hyb25vbG9naWN6bmllLg0KYGBge3J9DQpzb3VyY2VEYXRhRmlsZU5hbWUgPC0gInNsZWR6aWUuY3N2Ig0KYGBgDQoNCkRva29uYWpteSB3c3TEmXBuZWdvIHdjenl0YW5pYSBkYW55Y2ggLSBwaWVyd3N6eWNoIDEwMCB3aWVyc3p5IHcgY2VsdSBhdXRvbWF0eWN6bmVnbyBva3JlxZtsZW5pYSBrbGFzeSBvcmF6IHBvZGdsxIVkdSBkYW55Y2gNCmBgYHtyfQ0KaW5pdGlhbEhlcnJpbmdEYXRhIDwtIHJlYWQuY3N2KHNvdXJjZURhdGFGaWxlTmFtZSwgbnJvd3MgPSAxMDApDQpjbGFzc2VzIDwtIHNhcHBseShpbml0aWFsSGVycmluZ0RhdGEsIGNsYXNzKQ0KY2xhc3Nlcw0KYGBgDQpKYWsgd2lkYcSHLCB3acSZa3N6b8WbxIcgZGFueWNoIHpvc3RhxYJhIHpha3dhbGlmaWtvd2FuYSBqYWtvIGZhY3RvciBtaW1vLCDFvGUgeiBvcGlzdSB3eW5pa2EsIMW8ZSBzxIUgdG8gbGljemJ5LiBXY3p5dGFqbXkgemF0ZW0gY2HFgnkgcGxpayB3IG9kcG93aWVkbmljaCBrbGFzYWNoIC0gd2FydG/Fm2NpIGJyYWt1asSFY2UgcmVwcmV6ZW50dWplIGA/YCwgbmF0b21pYXN0IGpha28gc2VwYXJhdG9yIGR6aWVzacSZdG55IHByenlqbXVqZW15IGtyb3BrxJkgYC5gLiBQbGlrIG5pZSB6YXdpZXJhIHLDs3duaWXFvCBrb21lbnRhcnp5Lg0KU3p5YmthIGFuYWxpemEgcGxpa3UgcG96d29sacWCYSByw7N3bmllxbwgbmEgaWRlbnR5ZmlrYWNqxJksIMW8ZSB6YmnDs3Igc2vFgmFkYSBzacSZIHogb2tvxYJvIDUyLjYgdHlzIHdpZXJzenkgLSBpbmZvcm1hY2rEmSB0xIUgd3lrb3J6eXN0YW15IHcgY2VsdSBwcnp5xZtwaWVzemVuaWEgaW1wb3J0dS4NCmBgYHtyfQ0KY2xhc3NlcyA8LSBjKCJOVUxMIiwgcmVwKGMoIm51bWVyaWMiKSwgMTMpLCAiZmFjdG9yIiwgIm51bWVyaWMiKQ0KaGVycmluZ0RhdGEgPC0gcmVhZC5jc3Yoc291cmNlRGF0YUZpbGVOYW1lLCBjb21tZW50LmNoYXIgPSAiIiwgaGVhZGVyPVRSVUUsIGNvbENsYXNzZXMgPSBjbGFzc2VzLCBucm93cz01MjYwMCwgZGVjID0gIi4iLCBuYS5zdHJpbmdzID0gYygiPyIpKQ0KYGBgDQoNCg0KIyBPbcOzd2llbmllIGRhbnljaA0KU3ByYXdkem15IGN6eSByb3ptaWFyeSB6YmlvcnUgcG9rcnl3YWrEhSBzacSZIHogb2N6ZWtpd2FuaWFtaQ0KYGBge3J9DQpkaW0oaGVycmluZ0RhdGEpDQpgYGANCg0KV3lnbMSFZGEgbmEgdG8sIMW8ZSBkYW5lIHpvc3RhxYJ5IHBvcHJhd25pZSB3Y3p5dGFuZS4NClpvYmFjem15IGplIHdpxJljDQpgYGB7cn0NCnByZXR0eVRhYmxlKGhlYWQoaGVycmluZ0RhdGEsIG4gPSAxMDApKQ0KYGBgDQoNCg0KTmFzdMSZcG5pZSBwcnplcHJvd2Fkem15IGtyw7N0a2llIHBvZHN1bW93YW5pZQ0KYGBge3J9DQpzdW1tYXJ5KGhlcnJpbmdEYXRhKQ0KYGBgDQpKYWsgd2lkYcSHLCDFm2xlZHppZSBtYWrEhSBkxYJ1Z2/Fm8SHIHcgcHJ6ZWR6aWFsZSAxOS0zMi41IGNtLCB6IG1lZGlhbsSFIHVzdGFsb27EhSBkbGEgd2FydG/Fm2NpIDI1LjUgY20sIG5handpxJljZWogxZtsZWR6aSB3ecWCYXdpYSBzacSZIHcgc2llcnBuaXUsIHBhxbpkemllcm5pa3Ugb3JheiBsaXBjdS4gTW/FvG5hIHLDs3duaWXFvCBzdHdpZXJkemnEhywgxbxlIGxpY3piYSBicmFrdWrEhWN5Y2ggZGFueWNoIGRsYSBww7NsICpjZmluKiwgKmNoZWwxKiwgKmxjb3AqIG9yYXogKnNzdCogamVzdCBwb3LDs3dueXdhbG5hOyBkbGEgaW5ueWNoIGNvbHVtbiB3c3p5c3RraWUgZGFuZSBzxIUgb2tyZcWbbG9uZS4gDQpNb8W8ZW15IHLDs3duaWXFvCB6d3LDs2NpxIcgdXdhZ8SZIG5hIHdpZWxrb8WbxIcgcG/Fgm93w7N3IC0gaWNoIGxpY3piYSBtaWXFm2NpIHNpxJkgdyBwcnplZHppYWxlIDE0LjQgLSAxMDEuNiB0eXMgeiBtZWRpYW7EhSBuYSBwb3ppb21pZSA1NCB0eXMgb3JheiDFm3JlZG5pxIUgNTEuNSB0eXMuIFJvemvFgmFkIHRlbiBqZXN0IHdpxJljIHByYXdkb3BvZG9ibmllIHpibGnFvG9ueSBkbyBub3JtYWxuZWdvLCBjbyB6IHJlc3p0xIUgbW/FvGVteSBzcHJhd2R6acSHLg0KDQpgYGB7cn0NCiAgZ2dwbG90KGRhdGE9aGVycmluZ0RhdGEsIGFlcyh4PXRvdGFsbikpICsNCiAgZ2VvbV9kZW5zaXR5KCkgKw0KICBsYWJzKHggPSAiV2llbGtvxZvEhyBwb8WCb3d1IiwgeSA9ICJHxJlzdG/Fm8SHIikgKw0KICB0aGVtZV9idygpDQpgYGANCk5pZWtvbmllY3puaWUgLSBuaWUgcHJ6eXBvbWluYSBvbiDFvGFkbmVnbyB6IHBvcHVsYXJueWNoIHJvemvFgmFkw7N3Lg0KDQpNYW15IHphIHphZGFuaWUgendlcnlmaWtvd2HEhywgY3p5IHByemV6IG9zdGF0bmllIGtpbGthIGxhdCByb3ptaWFyIMWbbGVkemlhIHNwYWTFgi4gUHJ6eWRhdG55IG1vxbxlIG9rYXphxIcgc2llIHdpxJljIHd5a3JlcyBwcmV6ZW50dWrEhWN5IMWbcmVkbmkgcm96bWlhciBuYSByb2suIEplc3QgamVkbmFrIGplZGVuIHByb2JsZW0gLSB3IGRhbnljaCBtYW15IGplZHluaWUgbWllc2nEhWMgcG/Fgm93dSwgaW5mb3JtYWNqxJkgbyB0eW0sIMW8ZSBzxIUgdG8gZGFuZSBjaHJvbm9sb2dpY3puaWUgcG9zb3J0b3dhbmUgb3JheiBwb2Nob2R6xIUgeiBvc3RhdG5pY2ggNjAgbGF0LiBOYWxlxbx5IHdpxJljIHcgbWlhcsSZIG1vxbxsaXdvxZtjaSBkb2RhxIcgcm9rIGRvIGRhbnljaCBiYXp1asSFYyBuYSB0eWNoIGluZm9ybWFjamFjaC4NClNwcmF3ZHpteSBqZWRuYWsgY3p5IG5pZSBtYSB3xZtyw7NkIHR5Y2ggZGFueWNoIGJyYWt1asSFY3ljaCB3YXJ0b8WbY2kNCmBgYHtyfQ0Kc3VtKGlzLm5hKGhlcnJpbmdEYXRhJHhtb250aCkpDQpgYGANCmJyYWsgLSBtb8W8ZW15IHdpxJljIHByemVqxZvEhyBkYWxlai4NCkRhbmUgc8SFIGNocm9ub2xvZ2ljem5lLCB0YWvFvGUgcG93aW5uacWbbXkgYnnEhyB3IHN0YW5pZSB3eXduaW9za293YcSHIHJvayBuYSBwb2RzdGF3aWUgem1pYW55IG1pZXNpxIVjYSBuYSBtbmllanN6eQ0KDQpgYGB7cn0NCmV2YWx1YXRlWWVhckJ5TW9udGhDaGFuZ2UgPC0gZnVuY3Rpb24oZmlzaGluZykgew0KICBjdXJZZWFyIDwtIDANCiAgcmVjZW50TW9udGggPC0gMA0KICBsIDwtIGxhcHBseShmaXNoaW5nLCBmdW5jdGlvbih2KSB7DQogICAgdjwtIGFzLmludGVnZXIodikNCiAgICBpZiAodiA8IHJlY2VudE1vbnRoKQ0KICAgICAgY3VyWWVhciA8PC0gY3VyWWVhciArIDENCiAgICByZWNlbnRNb250aCA8PC0gdg0KICAgIGN1clllYXINCiAgfSkNCiAgdW5saXN0KGwpDQp9DQoNCmhlcnJpbmdEYXRhV2l0aFllYXJCeU1vbnRoQ2hhbmdlIDwtIG11dGF0ZShoZXJyaW5nRGF0YSwgeWVhcj1ldmFsdWF0ZVllYXJCeU1vbnRoQ2hhbmdlKHhtb250aCkpDQpoZXJyaW5nRGF0YVdpdGhZZWFyQnlNb250aENoYW5nZQ0KYGBgDQpQb2dydXB1am15IGplIHBvIHJva3UgZHppxJlraSBmdW5rY2ppDQpgYGB7cn0NCnN0YXRzSW5ZZWFyIDwtIGZ1bmN0aW9uKGRhdGFXaXRoWWVhcikgew0KICBkYXRhV2l0aFllYXIgJT4lDQogICAgZ3JvdXBfYnkoeWVhcikgJT4lDQogICAgc3VtbWFyaXNlX2F0KHZhcnMobGVuZ3RoKSwgZnVucyhtZWFuKC4sIG5hLnJtPVRSVUUpLCBtaW4gPSBtZWFuIC0gc2QoLiwgbmEucm09IFRSVUUpLCBtYXggPSBtZWFuICsgc2QoLiwgbmEucm09IFRSVUUpKSkNCn0NCg0KbWVhbkxlbmd0aEluWWVhckNoYXJ0IDwtIGZ1bmN0aW9uKGRhdGFXaXRoWWVhcikgew0KICBtZWFuUGVyWWVhciA8LSBzdGF0c0luWWVhcihkYXRhV2l0aFllYXIpDQogIGdncGxvdChkYXRhPW1lYW5QZXJZZWFyLCBhZXMoeCA9IHllYXIsIHkgPSBtZWFuKSkgKw0KICAgIGdlb21fcmliYm9uKGFlcyh5bWluID0gbWluLCB5bWF4ID0gbWF4KSwgYWxwaGEgPSAwLjUsIGZpbGwgPSAiZGFya3NlYWdyZWVuMyIsIGNvbG9yID0gInRyYW5zcGFyZW50IikgKw0KICAgICAgZ2VvbV9saW5lKGNvbG9yID0gImFxdWFtYXJpbmU0IiwgbHdkID0gMC43KSArIA0KICAgIGxhYnMoeCA9ICJSb2siLCB5ID0gIkTFgnVnb8WbxIcgxZtsZWR6aWEiKSArDQogIHRoZW1lX2J3KCkNCn0NCmBgYA0KDQpgYGB7cn0NCm1lYW5MZW5ndGhJblllYXJDaGFydChoZXJyaW5nRGF0YVdpdGhZZWFyQnlNb250aENoYW5nZSkNCmBgYA0KWmHFgm/FvGVuaWUgbyBjaHJvbm9sb2dpaSBqYWsgd2lkYcSHIG9rYXphxYJvIHNpxJkgYsWCxJlkbmUsIHByenluYWptbmllaiB3IHpha3Jlc2llIGtvbGVqbm/Fm2NpIG1pZXNpxJljeS4gSmVzdCByw7N3bmllxbwgcG9sZSBgcmVjcmAgbcOzd2nEhWNlIG8gcG/Fgm93aWUgdyByb2t1LiBXaWR6aW15LCB6ZSB3eXN0xJlwdWplIHBld25hIHBvd3RhcnphbG5vxZvEhywgd2nEmWMgamXFm2xpIGRhbmUgYnnFgnlieSBjaHJvbm9sb2dpY3puZSwgYSB6YXJhemVtIGljaCByb3pwacSZdG/Fm8SHIGplc3QgY2HFgmtpZW0gc3BvcmEsIHBvd2lubm8gdG8gcG96d29sacSHIG5hIHpncnVwb3dhbmllLg0KYGBge3J9DQpldmFsdWF0ZVllYXIgPC0gZnVuY3Rpb24oZmlzaGluZykgew0KICBjdXJZZWFyIDwtIDANCiAgcmVjZW50WWVhckZpc2hpbmcgPC0gMA0KICBsIDwtIGxhcHBseShmaXNoaW5nLCBmdW5jdGlvbih2KSB7DQogICAgaWYgKHYgIT0gcmVjZW50WWVhckZpc2hpbmcpDQogICAgICBjdXJZZWFyIDw8LSBjdXJZZWFyICsgMQ0KICAgIHJlY2VudFllYXJGaXNoaW5nIDw8LSB2DQogICAgY3VyWWVhcg0KICB9KQ0KICB1bmxpc3QobCkNCn0NCg0KaGVycmluZ0RhdGFXaXRoWWVhckJ5UmVjciA8LSBtdXRhdGUoaGVycmluZ0RhdGEsIHllYXI9ZXZhbHVhdGVZZWFyKHJlY3IpKQ0KaGVycmluZ0RhdGFXaXRoWWVhckJ5UmVjcg0KYGBgDQpgYGB7cn0NCm1lYW5MZW5ndGhJblllYXJDaGFydChoZXJyaW5nRGF0YVdpdGhZZWFyQnlSZWNyKQ0KYGBgDQpOaWVzdGV0eSwgcG9ub3duaWUgcG9yYcW8a2EuIFcgdGFraW0gcmF6aWUgbW/FvGUgcG8gcHJvc3R1IHNwcmF3ZHpteSBjenkgd3lzdMSZcHVqZSB0cmVuZCB3IGN6eXN0eWNoIGRhbnljaDoNCmBgYHtyfQ0KZ2dwbG90KGhlcnJpbmdEYXRhLCBhZXMoeD1hcy5udW1lcmljKHJvdy5uYW1lcyhoZXJyaW5nRGF0YSkpLCB5PWhlcnJpbmdEYXRhJGxlbmd0aCkpICsNCiAgZ2VvbV9saW5lKGNvbG9yID0gImFxdWFtYXJpbmU0IikgKw0KICBsYWJzKHggPSAiUm9rIiwgeSA9ICJExYJ1Z2/Fm8SHIMWbbGVkemlhIikgKyANCiAgdGhlbWVfYncoKQ0KYGBgDQpEb2JyemUsIGpha8SFxZsgemFsZcW8bm/Fm8SHIHdpZGHEhyAtIG1vxbxlIGFieSBiecSHIGplaiBwZXdueW0gdcWbcmVkbmlqbXkgZGFuZSBzdG9zdWrEhWMgcHJhd28gd2llbGtpY2ggbGljemIgb3JheiB3aWVkesSZIG8gdHltLCDFvGUgbWFteSBkbyBjenluaWVuaWEgeiBkYW55bWkgeiBva3Jlc3UgNjAgbGF0DQpgYGB7cn0NCmV4cGVjdGVkWWVhcnMgPC0gNjANCmV4cGVjdGVkUm93c1BlclllYXJzIDwtIG5yb3coaGVycmluZ0RhdGEpICUvJSA2MA0KaGVycmluZ0RhdGFXaXRoWWVhckJ5UGVyaW9kIDwtIGhlcnJpbmdEYXRhICU+JQ0KICBtdXRhdGUoeWVhciA9IDE6bigpICUvJSBleHBlY3RlZFJvd3NQZXJZZWFycykNCmhlcnJpbmdEYXRhV2l0aFllYXJCeVBlcmlvZA0KYGBgDQpgYGB7cn0NCm1lYW5MZW5ndGhJblllYXJDaGFydChoZXJyaW5nRGF0YVdpdGhZZWFyQnlQZXJpb2QpDQpgYGANCldpZHppbXkgd2nEmWMsIMW8ZSB3IHJ6ZWN6eXdpc3RvxZtjaSBpc3RuaWVqZSBpc3RvdG55IHNwYWRlayDFm3JlZG5pZWogZMWCdWdvxZtjaSDFm2xlZHppYSBuYSBwcnplxYJvbWllIGxhdC4gTW/FvGVteSB0byB6d2l6dWFsaXpvd2HEhy4NCmBgYHtyfQ0KaW1nIDwtIGltYWdlX3JlYWQoImhlcnJpbmcuc3ZnIikNCmltZyA8LSBpbWFnZV90cmltKGltZykNCmhlaWdodFJhdGlvIDwtIDExMC80NzUNCm1heFdpZHRoID0gNTQyDQptYXhIZWlnaHQgPC0gaGVpZ2h0UmF0aW8gKiBtYXhXaWR0aA0KaGVycmluZ3NJbWdzIDwtIHN0YXRzSW5ZZWFyKGhlcnJpbmdEYXRhV2l0aFllYXJCeVBlcmlvZCkgJT4lDQogIGFwcGx5KDEsIGZ1bmN0aW9uKHIpIHsNCiAgICBuZXdXaWR0aCA8LSByb3VuZCgocltbMl1dIC0gMTApICogMzApDQogICAgd2lkdGhEaWYgPC0gcm91bmQoKG1heFdpZHRoIC0gbmV3V2lkdGgpLzIpDQogICAgaGVpZ2h0RGlmIDwtIHJvdW5kKHdpZHRoRGlmICogaGVpZ2h0UmF0aW8pDQogICAgaW1hZ2Vfc2NhbGUoaW1nLCBuZXdXaWR0aCkgJT4lDQogICAgICAgIGltYWdlX2JvcmRlcigid2hpdGUiLCBwYXN0ZShjKHdpZHRoRGlmLCAieCIsIGhlaWdodERpZiksIGNvbGxhcHNlID0gIiIpKSAlPiUNCiAgICAgICAgaW1hZ2VfYW5ub3RhdGUocGFzdGUoYygiUm9rOiAiLCByW1sxXV0pLCBjb2xsYXBzZSA9ICIiKSwgY29sb3IgPSAiYmxhY2siLCBzaXplID0gMzApDQogIH0pDQpmcmFtZXMgPC0gaW1hZ2VfbW9ycGgoaW1hZ2Vfam9pbihoZXJyaW5nc0ltZ3MpLCBmcmFtZXMgPSAxKQ0KYW5pbSA8LSBpbWFnZV9hbmltYXRlKGZyYW1lcywgZnBzID0gMTApDQppbWFnZV93cml0ZShhbmltLCAiaGVycmluZy5naWYiKQ0KYGBgDQoNCiFbXShoZXJyaW5nLmdpZikNCg0KWmlkZW50eWZpa3VqbXkgbW/FvGxpd2UgcG93b2R5IHRha2llZ28gemFjaG93YW5pYSB3eXpuYWN6YWrEhWMgbWFjaWVyeiBrb3JlbGFjamkuDQoNCldjemXFm25pZWogamVkbmFrIG5hbGXFvHkgem1pZW5pxIcgdHlwIGtvbHVtbnkgYHhtb250aGAgeiBwb3dyb3RlbSBuYSBsaWN6YsSZIGkgcHJ6eWdvdG93YcSHIGRhbmUuDQpgYGB7cn0NCmhlcnJpbmdEYXRhV2l0aE51bWVyaWNNb250aCA9IG11dGF0ZShoZXJyaW5nRGF0YSwgeG1vbnRoID0gYXMubnVtZXJpYyh4bW9udGgpKQ0KYGBgDQpab2JhY3pteSB3eW5payB3IGJhcmR6aWVqIHByenlqYXpuZWogZm9ybWllDQpgYGB7cn0NCmdnY29ycigNCiAgaGVycmluZ0RhdGFXaXRoTnVtZXJpY01vbnRoLA0KICBuYnJlYWtzID0gOSwNCiAgbGFiZWwgPSBUUlVFLA0KICBsYWJlbF9zaXplID0gMywNCiAgY29sb3IgPSAiZ3JleTUwIiwNCiAgbGF5b3V0LmV4cCA9IDEsDQogIGhqdXN0ID0gMC43NSwNCikNCmBgYA0KV2lkYcSHIGR1xbzEhSBrb3JlbGFjasSZIHBvbWnEmWR6eSBwYXJhbWkgYGNoZWwxYCAtIGBsY29wMWAsIGBjaGVsMmAgLSBgbGNvcDJgLCBgZmJhcmAgLSBgY3VtZmAsIGBjZmluMmAgLSBgbGNvcDJgIG9yYXogYGN1bWZgIC0gYHRvdGFsbmAuIE5hIHNhbcSFIGTFgnVnb8WbxIcgdyBuYWp3acSZa3N6eW0gc3RvcG5pdSB3cMWCeXdhIGBzc3RgLCBgbmFvYCBvcmF6IGBmYmFyYC4NCldpZHppbXkgxbxlIGBsY29wMWAsIGBsY29wMmAgb3JheiBgY3VtZmAgbW9nxIUgem9zdGHEhyBwb21pbmnEmXRlIHcgcHJvY2VzaWUgZGFsc3plaiBhbmFsaXp5IGpha28gcmVkdW5kYW50bmUuIFRlb3JldHljem5pZSBtb2dsaWJ5xZtteSByw7N3bmllxbwgdXN1bsSFxIcgbWllc2nEhWMgemUgd3pnbMSZZHUgbmEgcG93dGFyemFsbm/Fm8SHLCBhbGUgbW/FvGUgbWllxIcgb24gd3DFgnl3IG5hIGN5a2wgcm96d29qb3d5LCB3acSZYyBwb3pvc3Rhd215IGdvIHDDs2tpIGNvLg0KYGBge3J9DQpmaWx0ZXJlZEhlcnJpbmdEYXRhV2l0aE51bWVyaWNNb250aCA8LSBzZWxlY3QoaGVycmluZ0RhdGFXaXRoTnVtZXJpY01vbnRoLCAtYyhsY29wMSwgbGNvcDIsIGN1bWYpKQ0KYGBgDQoNCg0KWm9iYWN6bXkgamVzemN6ZSBqYWsgcHJ6ZWJpZWdhIHphbGXFvG5vxZvEhyBwb21pxJlkenkgbmFqYmFyZHppZWogb2JpZWN1asSFY3ltaSB6bWllbm55bWkgbmEgZMWCdWdvxZvEhzsga29sb3JhbWkgd3lyw7PFvG5pb25vIHLDs3duaWXFvCBtaWVzacSFY2UNCmBgYHtyfQ0KbnVtZXJpY0NoYXJ0IDwtIGZ1bmN0aW9uKGNvbE5hbWUpIHsNCiAgZ2dwbG90KGhlcnJpbmdEYXRhLCBhZXMoeCA9IGxlbmd0aCwgeSA9IGhlcnJpbmdEYXRhW1tjb2xOYW1lXV0sIGNvbG9yPXhtb250aCkpICsNCiAgICBnZW9tX2ppdHRlcihzaXplID0gMS41LCBzdGF0ID0gImlkZW50aXR5IikgKw0KICAgIGxhYnMoeCA9ICJExYJ1Z2/Fm8SHIiwgeSA9IGNvbE5hbWUpDQp9DQpgYGANCg0KYGBge3J9DQpwbG90X2dyaWQoDQogIG51bWVyaWNDaGFydCgiY2hlbDEiKSwNCiAgbnVtZXJpY0NoYXJ0KCJmYmFyIiksDQogIG51bWVyaWNDaGFydCgic3N0IiksDQogIG51bWVyaWNDaGFydCgibmFvIikNCikNCmBgYA0KDQpQb3dpxIV6YW5pZSByb3NuxIVjZWdvIGBzc3RgIG9yYXogYG5hb2AgamFzbm8gd3DFgnl3YSBuZWdhdHl3bmllIG5hIGTFgnVnb8WbxIcsIHByenkgYGZiYXJgIHLDs3duaWXFvCBtb8W8bmEgZG9zdHJ6ZWMgcG96eXd5dG7EhSBrb3JlbGFjasSZLg0KDQpTcHJhd2R6bXkgcsOzd25pZcW8IHphbGXFvG5vxZvEhyBwb21pxJlkenkgbWllc2nEhWNlbSBwb8WCb3d1IGEgZMWCdWdvxZtjacSFOyB3IG1hY2llcnp5IGtvcmVsYWNqaSBuaWUgbW9nbGnFm215IHRlZ28gZG8ga29uY2EgendlcnlmaWtvd2HEhywgamVkbmFrIHplIHd6Z2zEmWR1IG5hIGN5a2xpY3pub8WbxIcgYnnFgm8gdG8gcsOzd25pZcW8IHV0cnVkbmlvbmUgKG1pZXNpxIVjIHcgdHltIHd5cGFka3UgamVzdCBzenR1Y3puaWUgbGljemLEhSkuDQoNCmBgYHtyfQ0KbnVtZXJpY0NoYXJ0KCJ4bW9udGgiKQ0KYGBgDQpXeWRhamUgc2nEmSwgxbxlIG5pZSBtYSDFvGFkbmVqIHpuYWN6xIVjZWogemFsZcW8bm/Fm2NpIHBvbWnEmWR6eSBtaWVzacSFY2VtIHBvxYJvd3UgYSBkxYJ1Z2/Fm2NpxIUgxZtsZWR6aWENCg0K